2011  Monitoring  Research  Review:  Ground-Based  Nuclear  Explosion  Monitoring  Technologies 


INFRASOUND  AS  A  DEPTH  DISCRIMINANT 

Stephen  J.  Arrowsmith,  Rod  W.  Whitaker,  and  Richard  J.  Stead 

Los  Alamos  National  Laboratory 
Sponsored  by  the  National  Nuclear  Security  Administration 

Award  No.  DE-AC52-06NA25396/LA09-Depth-NDD02 


ABSTRACT 


There  is  a  need  for  routine  measurements  of  infrasound  from  earthquakes,  in  conjunction  with  modeling,  to  better 
constrain  our  understanding  of  the  generation  of  infrasound  from  earthquakes,  in  particular  the  effect  of  source 
depth.  Here,  we  first  outline  a  systematic  search  for  infrasound  from  earthquakes  from  a  range  of  magnitudes.  Based 
on  our  observations,  earthquakes  with  Mw>~5  appear  to  be  reliable  sources  of  infrasound,  independent  of  depth, 
with  the  detection  of  infrasound  at  regional  distances  largely  driven  by  path  and  receiver  effects.  We  identified  clear 
infrasound  detections  from  large  earthquakes  as  deep  as  123  km.  In  contrast,  earthquakes  with  Mw<~4  are  poor 
infrasound  sources  and  are  unlikely  to  be  detected  regionally.  In  particular,  we  note  that  there  has  yet  to  be  a 
conclusive  detection  from  an  earthquake  with  Mw<3.0.  Our  results  imply  that  magnitude,  rather  than  depth,  appears 
to  be  the  dominant  factor  influencing  the  generation  of  infrasound  from  earthquakes.  To  understand  these 
observations  further,  we  have  developed  a  coupled  seismoacoustic  modeling  technique  comprising  a 
three-dimensional  (3D)  finite  difference  algorithm,  Rayleigh  integral,  and  time  domain  parabolic  equation.  By 
gluing  these  techniques  together,  we  are  modeling  the  generation  and  propagation  of  infrasound  from  earthquakes. 
The  effect  of  depth  on  these  simulations,  with  comparison  to  our  observational  dataset,  is  discussed. 
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OBJECTIVES 

The  goal  of  this  research  is  to  determine  whether  infrasound  might  be  useful  for  separating  shallow  events  from 
deep  events  and  possibly  providing  some  estimate  of  focal  depth.  Our  findings  address  this  question  but  also  provide 
more-general  insights  into  the  generation  of  infrasound  from  earthquakes.  By  modeling  the  generation  of  infrasound 
from  earthquakes  using  a  coupled  seismoacoustic  numerical  approach,  we  show  that  the  infrasound  amplitude 
should  indeed  scale  with  earthquake  depth,  with  a  linear  relation  between  the  log  of  the  infrasonic  amplitude  and  the 
log  of  the  depth.  However,  these  findings  are  not  borne  out  in  the  data  acquired  during  this  project,  which 
demonstrate  a  good  relationship  between  earthquake  magnitude  and  the  amplitude,  but  no  clear  relationship  with 
depth.  These  findings  highlight  the  need  for  more  infrasound  measurements  of  earthquakes  at  different  depths  in 
order  to  uniquely  separate  the  various  source  and  path  effects  in  the  real  data. 

RESEARCH  ACCOMPLISHED 


Data  Analysis 


The  development  of  an  empirical  dataset  is  well  described  in  the  previous  MRR  papers  under  this  project  and  is 
therefore  not  discussed  in  detail  here  (Arrowsmith  et  ah,  2010).  The  lack  of  adequate  data  has  hampered  this  effort, 
and  it  has  been  necessary  to  revise  our  strategy  from  requiring  multiple  array  detections  for  robust  association  to 
focusing  on  signals  from  large  earthquakes  recorded  within  1  skip  distance  in  order  to  maintain  an  adequate 
sampling.  The  final  dataset  is  summarized  in  Figure  1,  which  shows  that  the  ambient  noise  levels  appear  to  drive 
whether  or  not  we  observe  detections,  with  no  clear  correlation  between  the  detectability  and  the  focal  depth.  To 
explore  these  relationships  further.  Figure  2  shows  a  comparison  between  the  seismic  magnitude  and  predicted  path 
loss  for  the  five  undetected  earthquakes  with  root  mean  square  (RMS)  noise  levels  below  the  dashed  line  in  Figure  1 
and  the  detected  events.  This  figure  shows  that  the  undetected  earthquakes  associated  with  low  RMS  background 
noise  are  less  likely  to  be  detected  due  to  their  smaller  size  and  high  predicted  path  loss. 
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Figure  1.  Summary  of  infrasonic  detections  from  the  earthquake  dataset.  An  arbitrary  threshold  of  50  km  is 
used  to  distinguish  between  “shallow”  and  “deep”  earthquakes.  Note  that  for  detected  events,  both 
an  RMS  noise  and  peak-to-peak  signal  value  are  plotted;  for  undetected  events,  only  an  RMS  noise 
value  is  plotted. 
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Category 


Detected? 

O  no  (&  low  noise) 
O  yes 


Figure  2.  Comparison  between  undetected  events  with  low  noise  (Figure  1)  and  detected  events  in  terms 

of  the  seismic  magnitude  and  predicted  path  loss  from  parabolic  equation  modeling.  Four  out  of  five 
of  the  undetected  events  have  comparatively  low  magnitudes  and  high  path  loss — consistent  with 
the  lack  of  detection.  However,  one  anomalous  event  has  low  path  loss  and  a  relatively  high 
magnitude  of  5.7. 


For  each  detected  earthquake,  we  have  measured  a  variety  of  signal  parameters  using  the  beamformed, 
instrument-corrected  trace.  These  parameters  include  the  peak-to-peak  amplitude  in  Pa,  the  period  at  maximum 
amplitude,  the  group  velocity,  the  backazimuth  deviation,  and  the  cross-correlation  coefficient. 

The  dominant  frequency  content  of  most  of  the  infrasound  signals  is  —  1  Hz.  Spectra  showing  the  detailed  frequency 
content  for  four  selected  earthquakes  are  plotted  in  Figure  3.  The  spectra  clearly  show  that,  while  the  spectral 
content  is  indeed  centered  on  -  1  Hz,  the  signals  are  typically  in  the  range  of  -0.5-3  Hz  (although  this  varies 
somewhat  from  event  to  event).  As  shown  in  Figure  1,  peak-to-peak  signal  amplitudes  span  two  orders  of  magnitude 
from  -0.02  to  2  Pa.  Ambient  noise  levels  at  the  arrays  suggest  that  the  robust  detection  of  signals  below  -0.02  Pa  is 
likely  to  be  difficult  (Figure  1).  Azimuthal  deviations  have  a  standard  deviation  of  5.9°,  broadly  consistent  with 
magnitudes  of  backazimuth  deviations  predicted  over  multiple  years  from  ray  tracing  (Drob  et  ah,  2010). 

Modeling 

We  have  developed  a  coupled  numerical  source  model  for  modeling  the  generation  of  infrasound  from  earthquakes. 
The  coupled  model  fuses  a  three-dimensional  (3D)  finite  difference  time  domain  code  (ELAS3D — Larson  and 
Shultz,  1995)  with  a  Rayleigh  integral  code  that  we  have  developed.  Acceleration  time-histories  are  calculated  on  a 
2D  grid  centered  on  the  epicenter,  providing  the  initial  conditions  for  the  Rayleigh  integral  calculation.  As  discussed 
in  detail  in  Blackstock  (2000),  the  Rayleigh  integral  can  be  used  to  calculate  the  acoustic  pressure  at  some  distance 
from  a  baffled  piston  (i.e.,  all  acoustic  radiation  is  constrained  to  propagate  in  the  forward  direction)  as  an  integral  of 
contributions  from  accelerations  of  all  infinitesimal  area  elements,  dS,  of  the  piston: 
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p(x,y,z;t) 


niip(x',y'\t-R/c0 ) 

Poj  2ttR 


where  x’,  y’  are  the  coordinates  of  the  source  point  on  the  piston,  R  is  the  distance  to  the  field  point  of  interest,  pQ  is 
the  air  density,  and  C0  is  the  speed  of  sound  in  air. 


Even!  1  Event  2 


Figure  3.  Spectra  for  infrasound  signals  (red  curves)  and  corresponding  pre-event  noise  (blue  curves)  for  four 
different  events.  The  top  panel  represents  two  relatively  low  seismic  magnitude  deep  events,  while 
the  bottom  panel  represents  two  shallow  earthquakes  from  Tanzania  recorded  at  I32KE. 


The  Rayleigh  integral  is  an  approximation  to  the  Helmholtz-Kirchoff  integral  theorem  and  is  strictly  valid  in  the 
far-field,  beyond  the  Rayleigh  distance,  defined  in  general  by 


R 


0 


S' 

I’ 


where  S  is  the  surface  area  of  the  piston.  This  limitation  imposes  a  constraint  on  the  area  of  the  piston  that  can  be 
modeled,  since  the  assumption  of  a  constant  sound  speed,  C0  ,  becomes  unrealistic  at  large  distances  from  the  source 

due  to  the  variation  of  temperature  and  winds  with  height.  We  have  performed  a  series  of  synthetic  calculations 
using  different  piston  areas,  and  we  find  that  the  far-field  approximation  is  appropriate  at  5  km  elevation  from  a 
3x3  km2  region  (Figure  4).  We  choose  this  as  an  approximate  limit  on  the  size  of  the  region  that  can  be  modeled. 

The  radiation  pattern  can  be  characterized  by  an  amplitude  directivity  factor,  D  ,  which  is  defined  as  the  pressure  at 
any  angle  6  relative  to  that  (at  the  same  range  r)  at  d  =  0: 
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D(e) 


P{r,d) 
P{r, 0) ' 


The  radiation  pattern  from  a  circular  piston  of  radius  a  can  be  calculated  from  (Blackstock,  2000) 


D(0) 


2Jx{kasvad) 

kasind 


where  Jx  is  the  first-order  Bessel  function  of  the  first  kind,  k  is  the  wavenumber  (k  =  2lC  I A),  and  d  is  the  angle 
from  the  zenith.  The  number  of  nulls  and  secondary  maxima  in  the  radiation  pattern  is  determined  by  the  size  of  ka. 
As  the  piston  size  increases,  relative  to  the  wavelength,  the  number  of  nulls  increases. 


Acoustic  Pressure  from  a  3x3  km2  region 


Figure  4.  Acoustic  pressure  from  a  3  x  3  km2  region  at  different  heights  (1-20  km  in  increments  of  1  km) 
assuming  a  constant  acceleration  of  1  m/s2.  Above  5  km,  the  predicted  pressure  falls  off  as  1/r, 
indicating  that  the  far-field  assumption  is  valid  above  this  height. 

The  January  3,  201 1,  Circleville  earthquake,  which  is  discussed  in  further  detail  in  Arrowsmith  et  al.  (2011,  this 
Proceedings)  provides  a  good  case  study  to  explore  the  effect  of  earthquake  depth  on  the  recorded  signals. 

Additional  details  on  this  event  are  available  in  Arrowsmith  et  al.  (201 1).  Using  the  Harvard  centriod  moment  tensor 
(CMT)  solution,  one-dimensional  (ID)  velocity  model  for  southern  Utah,  and  a  Gaussian  source  time  function  of  1  s 
duration,  we  model  the  infrasound  generated  due  to  accelerations  in  the  3  x  3  km2  region  centered  on  the  epicenter 
and  calculate  the  associated  radiation  pattern  at  a  frequency  of  1  Hz  (Figure  5).  By  performing  ray-tracing  using  the 
G2S  atmospheric  model  for  this  event,  we  find  that  the  minimum  launch  angle  required  to  produce  a  return  for  the 
source-receiver  geometries  is  52°  from  the  zenith,  corresponding  to  a  maximum  value  of  D{9)  equal  to  0.005 

(Figure  5).  At  the  focal  depth  of  12  km,  we  estimate  a  source  amplitude  of  332.8  Pa  on-axis,  which  equates  to  a 
maximum  of  1.66  Pa  at  the  source  (332.8  x  0.005)  that  can  be  ducted.  By  scaling  our  synthetics  accordingly,  and 
using  a  time-domain  parabolic  equation  (PE)  method  to  estimate  the  resultant  acoustic  pressures  at  each  receiver,  we 
find  remarkably  good  agreement  (the  average  difference  between  observed  and  predicted  pressure  is  0.05  Pa). 
Having  obtained  confidence  in  this  modeling  approach,  we  explore  the  effect  of  depth  on  the  near-field  simulations, 
noting  that  these  results  are  specific  to  a  particular  geology  and  source  mechanism.  The  acoustic  waveforms  for  four 
different  simulations  at  depths  of  5  km,  12  km  (the  Harvard  CMT  depth),  50  km,  and  100  km  (with  all  other 
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parameters  held  constant)  are  shown  in  Figure  6.  The  modeling  predicts  that  there  should  be  a  clear  relationship 
between  earthquake  depth  and  the  near-field  infrasound  amplitude  (Figure  7).  Of  course,  the  actual  amplitudes 
detected  in  the  far-field  would  be  reduced  by  directivity  and  path  effects,  but  this  factor  would  be  constant  assuming 
the  same  time  and  source-receiver  geometry. 


Figure  5.  Radiation  pattern  for  the  piston  modeled  for  the  Circleville  earthquake  at  a  frequency  of  1  Hz 
(wavelength  of  340  m). 


Scaling  Laws 

The  development  of  an  empirical  relationship  between  event  depth  and  the  far-field  infrasound  amplitude  requires 
some  correction  for  the  effect  of  path  on  the  reduction  of  infrasound  amplitudes.  Here,  we  explore  two  techniques 
for  performing  path  corrections.  The  first  technique  is  based  on  an  empirical  formula  derived  from  observations  of 
underground  nuclear  tests  by  Mutschlecner  et  al.  (1999)  and  provides  an  amplitude,  An  ,  normalized  for  distance 
and  stratospheric  wind  speeds: 


4.  =V0 


-kVd 


where  s  =1.45,  k  =0.018  s/m,  R  is  the  great  circle  distance  to  the  epicenter,  Rs  is  an  arbitrary  standard  distance, 
and  V,  is  the  stratospheric  wind  component  directed  from  source  to  array. 


The  second  technique  uses  a  PE  model  (West  et  al.,  1992)  to  predict  path  losses  due  to  both  geometric  spreading  and 
absorption  from  source  to  receiver.  The  PE  model  can  predict  a  path-loss  in  dB  relative  to  1  km,  which  can 
subsequently  be  used  to  scale  all  recorded  amplitudes  to  1  km  from 
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Mutschlecner  et  al.  (1999)  did  not  have  access  to  modem  four-dimensional  atmospheric  models  such  as  the  ground- 
to-space  (G2S)  model;  thus,  we  intend  to  determine  whether  the  second  technique  might  reduce  the  variance  in 
earthquake  scaling  laws.  One  notable  difference  between  the  two  techniques  is  that  the  PE  method  can  be  applied  to 
any  phase,  whereas  the  Mutschlecner  et  al.  (1999)  method  is  strictly  for  stratospheric  returns. 


D**)th  ■  100  km 


Figure  6.  Simulated  on-axis  near-field  infrasound  waveforms  at  5  km  above  the  earthquake  epicenter  from  an 
earthquake  at  5  km  depth  (top  left),  12  km  depth  (bottom  left),  50  km  depth  (top  right),  and  100  km 
depth  (bottom  right). 
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Figure  7.  Relationship  between  focal  depth  and  the  peak  modeled  near-field  on-axis  infrasonic  amplitude  for 
the  January  3,  2011,  Circleville  earthquake  (the  actual  earthquake  focal  depth  was  ~12  km). 


Following  Mutschlecner  and  Whitaker  (2005)  and  Le  Pichon  et  al.  (2006),  we  expect  to  obtain  a  linear  relationship 
between  magnitude  and  the  log  of  the  amplitude  (corrected  for  path  effects).  Our  results  (Figure  8)  show  that  indeed 
there  is  such  a  relationship,  albeit  with  a  relatively  low  correlation.  The  correlation  reduces  somewhat  after 
correcting  for  path  effects  using  PE  modeling  (i.e.,  reducing  all  amplitudes  to  1  km),  but  there  is  still  some  large 
variability  that  remains  unaccounted  for.  Using  the  Mutschlecner  et  al.  (1999  )  path  correction,  and  excluding 
arrivals  with  celerities  outside  the  stratospheric  range  of  0.28-0.31  km/s,  we  obtain  a  much  better  fit  (R2  =  0.65), 
albeit  with  fewer  data  points.  The  regression  obtained  is  similar  to  the  Mutschlecner  and  Whitaker  (2005) 
regression,  which  used  ML  instead  of  Mw: 


log(A„)  =  0.70Mw-5.96. 


Our  initial  hypothesis  to  test  was  that  the  scatter  in  the  plots  in  Figure  8  could  be  explained  by  differences  in  event 
depth.  However,  plotting  the  residuals  from  the  regression  against  depth  (Figure  9)  shows  that  there  is  no  clear 
relationship  with  depth.  Clearly,  one  issue  is  that  the  seismically  obtained  depths  are  themselves  poorly  constrained. 
However,  one  can  think  of  the  events  as  being  robustly  separated  into  two  well-constrained  populations — deep  and 
shallow.  As  shown  in  Figure  9,  there  is  no  separation  between  the  two  populations. 

A  Brief  Discussion  of  Small  Earthquakes 

This  paper  has  focused  on  the  study  of  large  earthquakes  with  Mw>5.0.  A  complementary  study  discussed  in 
Arrowsmith  et  al.  (2011,  this  Proceedings)  focuses  on  small  earthquakes,  which  appear  to  be  relatively  poor 
infrasound  sources.  At  low  magnitudes  (Mw<4.0),  we  find  that — in  comparison  with  manmade  surface  explosions 
of  equivalent  seismic  magnitude — earthquakes  are  rarely  detected  at  regional  distances.  This  point  is  highlighted 
here  because  it  suggests  that  the  use  of  infrasound  as  a  source-type  discriminant  may  be  of  more  practical  value  in 
the  near-term  than  the  use  of  infrasound  as  a  depth  discriminant. 
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Figure  8.  Graphs  showing  the  relationship  between  the  raw  peak-to-peak  amplitude  and  Mw  (top  panel),  the 
peak-to-peak  amplitude  with  the  Mutschlecner  et  al.  (1999)  stratospheric  wind  correction,  for 
signals  with  stratospheric  celerities  only  (center  panel)  and  the  peak-to-peak  amplitude — corrected 
for  path  loss — and  Mw  (bottom  panel).  We  note  that  similar  plots  were  obtained  for  our 
measurements  of  RMS  amplitudes  that  were  evaluated  in  10  s  windows  centered  on  the  peak 
amplitude.  Black  triangles  depict  measurements  of  signals  from  the  2011  Circleville,  Utah, 
earthquake.  The  red  lines  on  each  plot  denote  the  best  least  squares  fit,  which  is  best  after  applying 
the  Mutschlecner  and  Whitaker  correction  (R2  =  0.65). 
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Figure  9.  Relationship  between  residuals  from  the  regression  shown  in  Figure  8  (lower  panel)  and  the  depth 
obtained  from  the  United  States  Geological  Survey  (USGS)  earthquake  catalog.  Although  seismic 
depths  may  be  poorly  constrained,  the  separation  between  shallow  (<  50  km)  and  deep  (>  50  km) 
should  be  robust.  Black  triangles  depict  measurements  of  signals  from  the  2011  Circleville,  Utah, 
earthquake. 


CONCLUSIONS  AND  RECOMMENDATIONS 


•  We  have  developed  a  numerical  approach  for  modeling  the  generation  of  infrasound  from  earthquakes  and 
show  that  for  the  January  3  Circleville,  Utah,  earthquake,  the  predicted  amplitudes  are  consistent  with  the 
observations. 

•  Numerical  modeling  of  the  generation  of  infrasound  from  earthquakes  predicts  a  linear  relationship 
between  the  log  of  the  infrasound  amplitude  and  the  log  of  depth. 

•  Empirical  analysis  confirms  that  infrasound  amplitude  scales  with  magnitude — as  previously  confirmed  by 
Mutschlecner  and  Whitaker  (2005)  and  Le  Pichon  et  al.  (2006) — but  indicates  that  depth  either  plays  no 
significant  role  or  else  is  obscured  by  other  effects  that  have  not  been  accounted  for  (in  contrast  to  the 
modeling). 

•  The  Mutschlecner  et  al.  (1999)  path  correction  appears  to  be  more  robust  for  predicting  path  corrections 
than  is  PE  modeling. 

•  Future  research  is  recommended  in  the  following  areas:  (1)  continued  dataset  development, 

(2)  improvements  to  source  modeling  (Kirchoff-Helmholtz),  and  (3)  improvements  to  path  corrections. 
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